Calibration Update
- Calibration period now from 2005-06-27 to 2015-01-01
- Spin up still 190. Should this stay the same?
- GLM not producing any ice at all. GOTM not following the same ice trend as other models.
- Flake is not below RMSE of 2, currently sitting at 2.33. It has very negative residuals around the thermocline depth, all the way out to -15. Are there any other parameters that can be changed for this model? If not, is this an adequate RMSE or does FLake have to be removed?
library(gotmtools)
## Loading required package: rLakeAnalyzer
## Warning: replacing previous import 'stats::filter' by 'dplyr::filter' when
## loading 'gotmtools'
library(LakeEnsemblR)
##
##
## _ _ _____ _ _ ____
## | | __ _| | _____| ____|_ __ ___ ___ _ __ ___ | |__ | | _ \
## | | / _` | |/ / _ | _| | '_ \/ __|/ _ | '_ ` _ \| '_ \| | |_) |
## | |__| (_| | | __| |___| | | \__ | __| | | | | | |_) | | _ <
## |_____\__,_|_|\_\___|_____|_| |_|___/\___|_| |_| |_|_.__/|_|_| \_\
##
##
## https://github.com/aemon-j/LakeEnsemblR
## LakeEnsemblR version 1.0.9 (2021-08-26) is loaded
## WARNING: Your current version of LakeEnsemblR (v1.0.9) is ahead of the master branch version (1.0.5)
## Development version may have an unexpected behaviour
##
## Attaching package: 'LakeEnsemblR'
## The following object is masked from 'package:gotmtools':
##
## analyse_strat
library(ggplot2)
library(LakeEnsemblR)
library(ggpubr)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(rLakeAnalyzer)
library(reshape)
##
## Attaching package: 'reshape'
## The following object is masked from 'package:dplyr':
##
## rename
library(reshape2)
##
## Attaching package: 'reshape2'
## The following objects are masked from 'package:reshape':
##
## colsplit, melt, recast
library(RColorBrewer)
library(lubridate)
##
## Attaching package: 'lubridate'
## The following object is masked from 'package:reshape':
##
## stamp
## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, union
library(Metrics)
##
## Attaching package: 'Metrics'
## The following object is masked from 'package:gotmtools':
##
## rmse
setwd("~/Dropbox/sunapee_LER_projections/LER_calibration/")
ncdf <- "output/ensemble_output_all_models_31Aug21.nc"
model <- c("FLake", "GLM", "GOTM", "Simstrat", "MyLake")
spin_up <- 190
plot_heatmap(ncdf, model = model) +
scale_colour_gradientn(limits = c(0, 32),
colours = rev(RColorBrewer::brewer.pal(11, "Spectral"))) + theme_classic()
## Extracted temp from output/ensemble_output_all_models_31Aug21.nc
## [,1]
## short_name "temp"
## units "degC"
## dimensions "lon lat z model member"
## Scale for 'colour' is already present. Adding another scale for 'colour',
## which will replace the existing scale.
plot_ensemble(ncdf, model = model, var = "ice_height")
fit <- calc_fit(ncdf, model = model, spin_up = spin_up)
fit
## $FLake
## rmse nse r bias mae nmae
## 1 2.324808 0.9060516 0.9571831 -0.7640633 1.477884 0.3071051
##
## $GLM
## rmse nse r bias mae nmae
## 1 1.714442 0.9470689 0.9736236 0.1523328 1.136328 0.3467074
##
## $GOTM
## rmse nse r bias mae nmae
## 1 1.85635 0.9379438 0.9700531 -0.3822369 1.26915 0.3418291
##
## $Simstrat
## rmse nse r bias mae nmae
## 1 1.50802 0.9590475 0.9812579 0.1904871 1.161026 0.2203282
##
## $MyLake
## rmse nse r bias mae nmae
## 1 1.409933 0.9642017 0.9820286 0.09948835 0.9535173 0.3288752
##
## $ensemble_mean
## rmse nse r bias mae nmae
## 1 1.17144 0.9752881 0.9878914 -0.1224119 0.7751371 0.2228305
# out <- analyze_ncdf(ncdf, model, spin_up = 190)
# out$stats
## Plot residuals
plist <- plot_resid(ncdf = "output/ensemble_output_all_models_14Aug21.nc", var = "temp")
## Extracted temp from output/ensemble_output_all_models_14Aug21.nc
## [,1]
## short_name "temp"
## units "degC"
## dimensions "lon lat z model member"
ggarrange(plotlist = plist)
Schmidt Stability
- Not much to say about Schmidt Stability, observations seem to track well with the mean at the very least
- RMSE shows that all models track relatively similarly in terms of GOF with the exception of FLake which has an excessively high RMSE at 239.7
- Flake Schmidt stability looks quite low after changing bathymetry to include the surface value only. If possible a double check on my code for this section would be much appreciated to ensure a mistake wasn’t made in the process.
ncdf <- "~/Dropbox/sunapee_LER_projections/LER_calibration/output/ensemble_output_all_models_31Aug21.nc"
######################## Calculating Schmidt Stability ################################
out <- load_var(ncdf = ncdf, var = "temp")
## [,1]
## short_name "temp"
## units "degC"
## dimensions "lon lat z model member"
# out <- as.data.frame(out[[1]])
bathy <- read.csv('~/Dropbox/sunapee_LER_projections/LER_inputs/sunapee_hypso.csv')
colnames(bathy) <- c("depths", "areas")
ts.sch <- lapply(out, function(x) {
ts.schmidt.stability(x, bathy = bathy, na.rm = TRUE)
})
## Reshape to data.frame
df <- melt(ts.sch, id.vars = 1)
colnames(df)[4] <- "model"
df$yday <- yday(df$datetime)
df$year <- year(df$datetime)
df <- filter(df, model != "FLake")
## CALCULATING FLAKE SCHMIDT STABILITY
listflake <- list(out[["FLake"]])
flake_depth <- glmtools::get_nml_value(nml_file = "~/Dropbox/sunapee_LER_projections/LER_calibration/FLake/flake.nml", arg_name = "depth_w_lk")
bathy <- filter(bathy, depths <= flake_depth)
ts.sch <- lapply(listflake, function(x) {
ts.schmidt.stability(x, bathy = bathy, na.rm = TRUE)
})
dfflake <- melt(ts.sch, id.vars = 1)
colnames(dfflake)[4] <- "model"
dfflake$model <- "FLake"
dfflake$yday <- yday(dfflake$datetime)
dfflake$year <- year(dfflake$datetime)
colnames(df)
## [1] "datetime" "variable" "value" "model" "yday" "year"
colnames(dfflake)
## [1] "datetime" "variable" "value" "model" "yday" "year"
df <- rbind(df, dfflake)
########################## END CALCULATING FLAKE SCHMIDT STABILITY
wideform <- dcast(df, datetime~model, value.var = "value")
wideform <- filter(wideform, is.na(Obs) == FALSE & is.na(GLM) == FALSE &
is.na(GOTM) == FALSE & is.na(FLake) == FALSE &
is.na(Simstrat) == FALSE & is.na(MyLake) == FALSE)
rmse <- c()
models <- c("FLake", "GOTM", "Simstrat", "MyLake", "GLM")
rmse <- c(rmse, (rmse(wideform$Obs, wideform$FLake)))
rmse <- c(rmse, (rmse(wideform$Obs, wideform$GOTM)))
rmse <- c(rmse, (rmse(wideform$Obs, wideform$Simstrat)))
rmse <- c(rmse, (rmse(wideform$Obs, wideform$MyLake)))
rmse <- c(rmse, (rmse(wideform$Obs, wideform$GLM)))
models <- as.data.frame(models)
rmse <- as.data.frame(rmse)
model_rmse <- cbind(models, rmse)
model_rmse
## models rmse
## 1 FLake 166.16345
## 2 GOTM 69.78399
## 3 Simstrat 75.21572
## 4 MyLake 55.53663
## 5 GLM 113.50085
dfobs <- filter(df, model == "Obs")
df <- filter(df, model != "Obs")
df <- df %>%
dplyr :: group_by(yday, year) %>%
dplyr :: mutate(mean = mean(value, na.rm = TRUE)) %>%
dplyr :: mutate(sd = sd(value, na.rm = TRUE))
dfobs$mean <- NA
dfobs$sd <- NA
df <- rbind(df, dfobs)
## plot results
ggplot(df, aes(yday, value, colour = model)) +
# geom_line(data=df, aes(y=mean, x=yday), color = "black", size = 1) +
geom_line() +
# geom_ribbon(data = df, aes(ymin = mean-sd, ymax = mean+sd), alpha = 0.6,
# linetype = 0.1,
# color = "grey") +
facet_wrap(~year) +
labs(y = "Schmidt stability (J/m2)") +
theme_classic() + ylim(-50, 1000)
ggplot(subset(df, model != "Obs"), aes(yday, mean)) +
geom_line() +
geom_ribbon(data = subset(df, model != "Obs"), aes(ymin = mean-sd, ymax = mean+sd), alpha = 0.6,
linetype = 0.1,
color = "grey") +
facet_wrap(~year) +
labs(y = "Schmidt stability (J/m2)") +
theme_classic() + ylim(-50, 1000) +
geom_line(data = subset(df, model == "Obs"), aes(yday, value, col = "Obs"))
Thermocline Depth
- Similar update on thermocline depth.
- GLM appears to be the poorest performer most years based on the first figure.
- The extra years of observations give me more confidence, as 2009 and 2010 were a bit off. The thermocline seems to fit quite well from 2011 and on.
- FLake is the highest performer for RMSE, with MyLake being the lowest performer. When looking at MyLake on the plots however, it seems to perform relatively well with occassional anomalous jumps down to the max depth towards the end of the season. GLM may perform worse than MyLake during the overall season as the RMSE is 7 and GLM consistently models significantly lower than the observed thermocline depth.
library(gotmtools)
library(LakeEnsemblR)
library(ggplot2)
library(LakeEnsemblR)
library(ggpubr)
library(dplyr)
library(rLakeAnalyzer)
library(reshape)
library(RColorBrewer)
library(lubridate)
ncdf <- "~/Dropbox/sunapee_LER_projections/LER_calibration/output/ensemble_output_all_models_31Aug21.nc"
out <- load_var(ncdf = ncdf, var = "temp")
## [,1]
## short_name "temp"
## units "degC"
## dimensions "lon lat z model member"
## Same for thermocline depth
ts.td <- lapply(out, function(x) {
ts.thermo.depth(x, Smin = 0.1, na.rm = TRUE)
})
df <- melt(ts.td, id.vars = 1)
colnames(df)[4] <- "model"
df$yday <- yday(df$datetime)
df$month <- month(df$datetime)
df$year <- year(df$datetime)
wideform <- dcast(df, datetime~model, value.var = "value")
wideform <- filter(wideform, is.na(Obs) == FALSE & is.na(GLM) == FALSE &
is.na(GOTM) == FALSE & is.na(FLake) == FALSE &
is.na(Simstrat) == FALSE & is.na(MyLake) == FALSE)
rmse <- c()
models <- c("FLake", "GOTM", "Simstrat", "MyLake", "GLM")
rmse <- c(rmse, (rmse(wideform$Obs, wideform$FLake)))
rmse <- c(rmse, (rmse(wideform$Obs, wideform$GOTM)))
rmse <- c(rmse, (rmse(wideform$Obs, wideform$Simstrat)))
rmse <- c(rmse, (rmse(wideform$Obs, wideform$MyLake)))
rmse <- c(rmse, (rmse(wideform$Obs, wideform$GLM)))
models <- as.data.frame(models)
rmse <- as.data.frame(rmse)
model_rmse <- cbind(models, rmse)
model_rmse
## models rmse
## 1 FLake 2.283110
## 2 GOTM 4.573914
## 3 Simstrat 3.415918
## 4 MyLake 13.902624
## 5 GLM 6.597361
dfobs <- filter(df, model == "Obs")
df <- filter(df, model != "Obs")
df <- df %>%
dplyr :: group_by(yday, year) %>%
dplyr :: mutate(mean = mean(value, na.rm = TRUE)) %>%
dplyr :: mutate(sd = sd(value, na.rm = TRUE))
dfobs$mean <- NA
dfobs$sd <- NA
df <- rbind(df, dfobs)
ggplot(subset(df, month >= 6 & month <= 8), aes(yday, value, colour = model)) +
facet_wrap(~year) +
geom_line() +
labs(y = "Thermocline depth (m)") +
scale_y_continuous(trans = "reverse") +
theme_classic()
ggplot(subset(df, month >= 6 & month <= 8 & model != "Obs"), aes(yday, mean)) +
facet_wrap(~year) +
geom_line() +
geom_ribbon(data = subset(df, month >= 6 & month <= 8 & model != "Obs"), aes(ymin = mean-sd, ymax = mean+sd), alpha = 0.6,
linetype = 0.1,
color = "grey") +
labs(y = "Thermocline depth (m)") +
scale_y_continuous(trans = "reverse") +
theme_classic() +
geom_line(data = subset(df, month >= 6 & month <= 8 & model == "Obs"), aes(yday, value, col = "Obs"))
dfsummer <- filter(df, month >=6 & month <=8)
All other metrics
- Observations are added for this group of metrics
- There seems to be an issue with 2009 calculations, as shown by the printed dataframe.
- The only metric to be compared with ice will have to be added in - coming soon.
- There is very limited stratification data calculated for the observed data.
- Seems to very high uncertainty for bottom temperature. Reason to remove it from desired metrics?
## [,1]
## short_name "temp"
## units "degC"
## dimensions "lon lat z model member"
## [,1]
## short_name "temp"
## units "degC"
## dimensions "lon lat z model member"
## [,1]
## short_name "ice_height"
## units "m"
## dimensions "lon lat model member"
## # A tibble: 24 x 6
## # Groups: year, variable [24]
## year variable value model mean sd
## <dbl> <fct> <dbl> <chr> <dbl> <dbl>
## 1 2009 TsMax 0.140 Obs NA NA
## 2 2009 TsMaxDay 5 Obs NA NA
## 3 2009 TsMin -0.270 Obs NA NA
## 4 2009 TsMinDay 0 Obs NA NA
## 5 2009 TbMax 0.970 Obs NA NA
## 6 2009 TbMaxDay 11 Obs NA NA
## 7 2009 TbMin 0.400 Obs NA NA
## 8 2009 TbMinDay 0 Obs NA NA
## 9 2009 MaxStratDur NA Obs NA NA
## 10 2009 MeanStratDur NA Obs NA NA
## # … with 14 more rows
Ice off
- GLM does not model ice so is not present
- FLake, MyLake and Simstrat do extremely well when comparing by plot and by RMSE, with all models hovering slightly above or below 5 RMSE.
- GOTM is performing poorly when modeling ice, with an RMSE of 59 and a consistently lower ice-off date than any of the other models.
## [,1]
## short_name "temp"
## units "degC"
## dimensions "lon lat z model member"
## [,1]
## short_name "temp"
## units "degC"
## dimensions "lon lat z model member"
## [,1]
## short_name "ice_height"
## units "m"
## dimensions "lon lat model member"
## models rmse
## 1 FLake 4.459696
## 2 GOTM 59.897134
## 3 Simstrat 4.333333
## 4 MyLake 4.887626
library(gotmtools)
library(LakeEnsemblR)
library(ggplot2)
library(LakeEnsemblR)
library(ggpubr)
library(dplyr)
library(rLakeAnalyzer)
library(reshape)
library(RColorBrewer)
library(lubridate)
setwd("~/Dropbox/sunapee_LER_projections/LER_calibration/calibration_results_FLake_083021_v2/")
out_f <- "calibration_results_FLake_083021_v2"
config_file <- 'LakeEnsemblRsun.yaml'
model <- c("FLake")
param_file <- "calibration_results_FLake_083021_v2/FLake_LHC_202108311031.csv"
cal_files <- list.files(".", full.names = TRUE)
cal_files <- cal_files[c(1,2)]
res <- load_LHC_results(config_file = config_file, model = model, res_files = cal_files)
dim(res[[model]])
## [1] 500 12
df <- plyr::ldply(res, function(x) {
df <- x[, -c(3:7)]
reshape2::melt(df, id.vars = c("par_id", "rmse"))
}, .id = "model")
df$id_no <- as.numeric(gsub(".*?([0-9]+).*", "\\1", df$par_id))
bst_par <- df$id_no[which.min(df$rmse)]
sub <- df[df$id_no == bst_par, ]
sub
## model par_id rmse variable value id_no
## 330 FLake p0330 2.324808 wind_speed 1.0404000 330
## 830 FLake p0330 2.324808 swr 0.8437100 330
## 1330 FLake p0330 2.324808 LAKE_PARAMS/c_relax_C 0.0054947 330
## 1830 FLake p0330 2.324808 LAKE_PARAMS/depth_bs_lk 2.4575000 330
## 2330 FLake p0330 2.324808 LAKE_PARAMS/T_bs_lk 15.2050000 330
p1 <- ggplot(df) +
geom_point(aes(value, rmse)) +
facet_wrap(model~variable, scales = "free_x") +
geom_hline(yintercept = 2, linetype = "dashed") +
ylab("RMSE (\u00B0C)") +
geom_vline(data = sub, aes(xintercept = value)) +
geom_hline(yintercept = 3.5, color = "red", linetype = "dashed") +
# coord_cartesian(ylim = c(1, 4)) +
# scale_x_log10() +
theme_classic(base_size = 16)
p1
library(gotmtools)
library(LakeEnsemblR)
library(ggplot2)
library(LakeEnsemblR)
library(ggpubr)
library(dplyr)
library(rLakeAnalyzer)
library(reshape)
library(RColorBrewer)
library(lubridate)
setwd("~/Dropbox/sunapee_LER_projections/LER_calibration/calibration_results_GLM_082621/")
out_f <- "calibration_results_GLM_082621"
config_file <- 'LakeEnsemblRsun.yaml'
model <- c("GLM")
# param_file <- "calibration_results_FLake_083021_v2/FLake_LHC_202108311031.csv"
cal_files <- list.files(".", full.names = TRUE)
cal_files <- cal_files[c(1,2)]
res <- load_LHC_results(config_file = config_file, model = model, res_files = cal_files)
dim(res[[model]])
## [1] 500 12
df <- plyr::ldply(res, function(x) {
df <- x[, -c(3:7)]
reshape2::melt(df, id.vars = c("par_id", "rmse"))
}, .id = "model")
df$id_no <- as.numeric(gsub(".*?([0-9]+).*", "\\1", df$par_id))
bst_par <- df$id_no[which.min(df$rmse)]
sub <- df[df$id_no == bst_par, ]
sub
## model par_id rmse variable value id_no
## 36 GLM p0036 1.714442 wind_speed 1.49590 36
## 536 GLM p0036 1.714442 swr 1.06550 36
## 1036 GLM p0036 1.714442 sediment/sed_temp_mean...3 3.48260 36
## 1536 GLM p0036 1.714442 sediment/sed_temp_mean...4 10.59700 36
## 2036 GLM p0036 1.714442 glm_setup/max_layer_thick 0.64125 36
p1 <- ggplot(df) +
geom_point(aes(value, rmse)) +
facet_wrap(model~variable, scales = "free_x") +
geom_hline(yintercept = 2, linetype = "dashed") +
ylab("RMSE (\u00B0C)") +
geom_vline(data = sub, aes(xintercept = value)) +
geom_hline(yintercept = 3.5, color = "red", linetype = "dashed") +
# coord_cartesian(ylim = c(1, 4)) +
# scale_x_log10() +
theme_classic(base_size = 16)
p1
## Warning: Removed 5 rows containing missing values (geom_point).
library(gotmtools)
library(LakeEnsemblR)
library(ggplot2)
library(LakeEnsemblR)
library(ggpubr)
library(dplyr)
library(rLakeAnalyzer)
library(reshape)
library(RColorBrewer)
library(lubridate)
setwd("~/Dropbox/sunapee_LER_projections/LER_calibration/calibration_results_GOTM_081321/")
out_f <- "calibration_results_GOTM_081321"
config_file <- 'LakeEnsemblRsun.yaml'
model <- c("GOTM")
# param_file <- "calibration_results_FLake_083021_v2/FLake_LHC_202108311031.csv"
cal_files <- list.files(".", full.names = TRUE)
cal_files <- cal_files[c(1,2)]
res <- load_LHC_results(config_file = config_file, model = model, res_files = cal_files)
dim(res[[model]])
## [1] 500 10
df <- plyr::ldply(res, function(x) {
df <- x[, -c(3:7)]
reshape2::melt(df, id.vars = c("par_id", "rmse"))
}, .id = "model")
df$id_no <- as.numeric(gsub(".*?([0-9]+).*", "\\1", df$par_id))
bst_par <- df$id_no[which.min(df$rmse)]
sub <- df[df$id_no == bst_par, ]
sub
## model par_id rmse variable value id_no
## 56 GOTM p0056 1.85635 wind_speed 5.4340e-01 56
## 556 GOTM p0056 1.85635 swr 1.3701e+00 56
## 1056 GOTM p0056 1.85635 turbulence/turb_param/k_min 2.6656e-06 56
p1 <- ggplot(df) +
geom_point(aes(value, rmse)) +
facet_wrap(model~variable, scales = "free_x") +
geom_hline(yintercept = 2, linetype = "dashed") +
ylab("RMSE (\u00B0C)") +
geom_vline(data = sub, aes(xintercept = value)) +
geom_hline(yintercept = 3.5, color = "red", linetype = "dashed") +
# coord_cartesian(ylim = c(1, 4)) +
# scale_x_log10() +
theme_classic(base_size = 16)
p1
library(gotmtools)
library(LakeEnsemblR)
library(ggplot2)
library(LakeEnsemblR)
library(ggpubr)
library(dplyr)
library(rLakeAnalyzer)
library(reshape)
library(RColorBrewer)
library(lubridate)
setwd("~/Dropbox/sunapee_LER_projections/LER_calibration/calibration_results_MyLake_081321/")
out_f <- "calibration_results_MyLake_081321"
config_file <- 'LakeEnsemblRsun.yaml'
model <- c("MyLake")
# param_file <- "calibration_results_FLake_083021_v2/FLake_LHC_202108311031.csv"
cal_files <- list.files(".", full.names = TRUE)
cal_files <- cal_files[c(1,2)]
res <- load_LHC_results(config_file = config_file, model = model, res_files = cal_files)
dim(res[[model]])
## [1] 500 10
df <- plyr::ldply(res, function(x) {
df <- x[, -c(3:7)]
reshape2::melt(df, id.vars = c("par_id", "rmse"))
}, .id = "model")
df$id_no <- as.numeric(gsub(".*?([0-9]+).*", "\\1", df$par_id))
bst_par <- df$id_no[which.min(df$rmse)]
sub <- df[df$id_no == bst_par, ]
sub
## model par_id rmse variable value id_no
## 219 MyLake p0219 1.409933 wind_speed 1.48920 219
## 719 MyLake p0219 1.409933 swr 1.27410 219
## 1219 MyLake p0219 1.409933 Phys.par/C_shelter 0.39897 219
p1 <- ggplot(df) +
geom_point(aes(value, rmse)) +
facet_wrap(model~variable, scales = "free_x") +
geom_hline(yintercept = 2, linetype = "dashed") +
ylab("RMSE (\u00B0C)") +
geom_vline(data = sub, aes(xintercept = value)) +
geom_hline(yintercept = 3.5, color = "red", linetype = "dashed") +
# coord_cartesian(ylim = c(1, 4)) +
# scale_x_log10() +
theme_classic(base_size = 16)
p1
library(gotmtools)
library(LakeEnsemblR)
library(ggplot2)
library(LakeEnsemblR)
library(ggpubr)
library(dplyr)
library(rLakeAnalyzer)
library(reshape)
library(RColorBrewer)
library(lubridate)
setwd("~/Dropbox/sunapee_LER_projections/LER_calibration/calibration_results_Simstrat_081321/")
out_f <- "calibration_results_Simstrat_081321"
config_file <- 'LakeEnsemblRsun.yaml'
model <- c("Simstrat")
# param_file <- "calibration_results_FLake_083021_v2/FLake_LHC_202108311031.csv"
cal_files <- list.files(".", full.names = TRUE)
cal_files <- cal_files[c(1,2)]
res <- load_LHC_results(config_file = config_file, model = model, res_files = cal_files)
dim(res[[model]])
## [1] 500 10
df <- plyr::ldply(res, function(x) {
df <- x[, -c(3:7)]
reshape2::melt(df, id.vars = c("par_id", "rmse"))
}, .id = "model")
df$id_no <- as.numeric(gsub(".*?([0-9]+).*", "\\1", df$par_id))
bst_par <- df$id_no[which.min(df$rmse)]
sub <- df[df$id_no == bst_par, ]
sub
## model par_id rmse variable value id_no
## 20 Simstrat p0020 1.502343 wind_speed 1.4883000 20
## 520 Simstrat p0020 1.502343 swr 1.0299000 20
## 1020 Simstrat p0020 1.502343 ModelParameters/a_seiche 0.0035929 20
p1 <- ggplot(df) +
geom_point(aes(value, rmse)) +
facet_wrap(model~variable, scales = "free_x") +
geom_hline(yintercept = 2, linetype = "dashed") +
ylab("RMSE (\u00B0C)") +
geom_vline(data = sub, aes(xintercept = value)) +
geom_hline(yintercept = 3.5, color = "red", linetype = "dashed") +
# coord_cartesian(ylim = c(1, 4)) +
# scale_x_log10() +
theme_classic(base_size = 16)
p1
